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METHOD OP SmPLIFYING THE TASK OP KODZIXXKG A POROUS 
GEOLOGICAL MEDltW CROSSED BY AK IRREGULAR 
FRACTURE NETWORK 

The invention relates to a method designed to 
simplify the task of modelling a porous geological 
medium- crossed by an irregular network of fissures, 
which makes it easier to correlate the characterisation 
of fractured reservoirs with dual-porosity models. The 
method may be used by reservoir engineers in oil 
production, for example, in order to produce reliable 
flow predictions. 

Fractured reservoirs are an extreme type of 
heterogeneous reservoirs made up of two contrasting 
media, a matrix medium containing the majority of the 
oil and of a low permeability and a fractured medium 
accounting for less than 1% of the oil in place and 
highly conductive. The fractured medium itself may be 
very complex having different fracture sets, each 
characterised by its density, length, orientation, tilt 
and aperture. 3D images of fractured reservoirs can not 
be used directly in the form of input data for 
simulating a reservoir. It has long been considered 
unrealistic to represent fracture networks within 
reservoir flow simulators because the layout of the 
network is partially unknown and because of the 
numerical limitations inherent in juxtaposing numerous 
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cells of extremely contrasting sizes and properties. A 
simplified hut realistic way of modelling such media 
would therefore be of great interest to reservoir 
engineers . 

The "dual porosity approach", as taught, for 
example, by Warren J.E, et al in "The Behaviour of 
Naturally Fractured Reservoirs", SPE Journal (September 
1963), 245-255, is known in the art as a means of 
interpreting single-phase flow behaviour observed when 
testing a fractured reservoir. In accordance with this 
basic model, any elementary volume of the fractured 
reservoir is modelled in the form of a set of 'identical 
parallelepipedic blocks limited by an orthogonal system 
of continuous uniform fractures oriented in the 
direction of one of the three main flow directions. On 
the reservoir scale, fluid flows through the fractured 
medium only and fluid exchanges occur locally between 
the fissures and the matrix blocks. 

Numerous simulators have been developed for 
fractured reservoirs using a model of this type with 
specific improvements in respect of the way flows are 
exchanged between matrix and fracture, governed by 
capillary, gravitational and viscous forces and 
compositional mechanisms as well as matrix-to-matrix 
flow exchanges {dual -permeability, dual-porosity 
simulators) . Several exati^)les of prior art techniques 
are cited in the reference works listed below. 
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Thomas, L.K. et al: -Fractured Reservoir 
Simulation", SPB Journal (February 1983) 42-54; 

Quandalle, P. et al : "Typical Features of a New 
Multipurpose Reservoir Simulator", SPE 16007, presented 
at the 9th Symposium on reservoir simulation held in 
San Antonio, Texas, 1-4 February" 1987; 

Coats, K. H. : "Implicit Compositional Simulation 
of Single-Porosity and Dual-Porosity Reservoirs", sPE 
18427 presented to the SPE Symposium on reservoir 
simulation held in Houston, Texas, 6-8 February 1989. 

One of the problems which reservoir engineers 
encounter is that of applying parameters to this basic 
model in order to produce reliable flow predictions. 
The permeabilities of equivalent fractures and the size 
of the matrix blocks in particular must be ascertained 
for each cell of the flow simulator. Whilst the 
permeability of the matrix can be estimated from cores, 
there is no simple method of estimating the 
permeability of the fracture network contained within 
the cell, i.e. the equivalent fracture permeability, 
and it is necessary to take account of the geometry and 
properties of the real fracture network. One method of 
determining equivalent fracture permeabilities in a 
fractured network is disclosed in patent application 
EN. 96/16330 filed in parallel with this one. 

A reference procedure is known for determining the 
dimensions a, b of each block of a section crossed by a 



18006661233; Mar^-00 9:47; 



regular grid of fractures Feq, which is equivalent to 
the section of a natural, fractured multi- layered 
medium crossed by a fracture network FN along a 
reference plane parallel with the layers (commonly a 
horizontal or essentially horizontal plane) . For each 
layer of the fractured rock volume studied (figure i) , 
the "horizontal" dimensions a, b of the blocks of the 
equivalent section are determined iteratively by means 
of computations and by comparing the oil recovery 
functions as a function of time R(t) and Req(t) 
respectively in the real section RE of the fractured 
rock volume studied and in the section EQ cut into 
identically sized "sugar lumps" equivalent to the 
distribution of real blocks. This conventional method 
requires a single -porosity multi-phase flow simulator 
which divides up the matrix blocks and fractures in 
such a way that the recovery curves can be compared. 
Such a procedure is extremely costly in view of the 
fact that the task of dividing up the real section may 
involve a very high number of cells. The real shape of 
the blocks in effect have to be represented using thin 
fracture cells along the boundaries of each block. The 
matrix also has to be divided up into a large enough 
number of cells to obtain an accurate block-fracture 
imbibition transfer function. 

Various prior art techniques used in the field are 
described, for example, in: 
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- Bourbiaux, B. et al : " Rxperimental Study of Cocurrent 
and Countercurrent Flows in Natural Porous Media", SPE 
Reservoir Engineering (August 1990), 361-368, 

- Cuiec, L. et al: "Oil Recovery by Imbibition in Low- 
Permeability Chalk", SPE Formation Evaluation 
(September 1994), 200-208. 

However, specific imbibition characteristics have 
never been used as a means of finding the dimensions of 
the equivalent block in dual porosity models. Reservoir 
engineers do not therefore have a systematic tool which 
can be used to calculate the dimensions of 
parallelepipedic blocks equivalent to mult i -phase flows 
in order to determine the real distribution of the 
blocks in each fractured reservoir zone. 
15 Techniques for integrating natural fracture data 

in fractured reservoir models are also known in the 
prior art. The fracture data are specifically data of a 
geometric nature and do not incorporate measurements of 
the density, length, azimuth and tilt of fracture 
planes, either as observed on outcrops, mine drifts and 
cores or inferred from well logging. Different fracture 
sets can be distinguished and characterised by 
different statistical distributions of their fracture 
attributes: Once the fracture patterns have been 
25 characterised, numerical networks of those fracture 
sets can be generated using a stochastic method which 
respects the statistical distributions of fracture 



20 
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parameters. Such methods are described in patents FR-A- 
2.725.814, 2.725.794 or 2.733.073, for example, filed 
by the applicant. 

The method of the invention provides a simplified 
5 approach to modelling an original geological medium 
which is heterogeneous and porous {such as a reservoir 
crossed by an irregular fracture network, for example) 
in the form of a transposed or equivalent medium in 
such a way that the transposed medium is equivalent to 
10 the original medium as regards a given physical 
transfer function which is known in the transposed 
medium. The method consists in: 

- forming an image of the geological medium in at least 
two dimensions in the form of an array of pixels and 

15 associating a specific initial value for said function 
to each pixel in the array, 

- determining step by step the value to be assigned to 
each pixel in the array for the physical transfer 
function, by reference to the function values assigned 

20 to adjacent pixels in the image, and 

- determining a physical property of the transposed or 
equivalent medium by identifying values of the transfer 
function known for the (simplified) transposed medium 
with the value deteirmined step by step for the original 

25 medium. 

The physical transfer function may represent 
variations between different parts of the geological 
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medium, variations in distances^ trauismiseivities or 
heat for example (such as transfers of heat between a 
reservoir and a well crossing through this reservoir) 
etc. 

The method can be applied as a means of 
determining from an image of a porous geological medium 
crossed by an irregular fracture network, for example, 
a transposed medium comprising a set of regularly- 
arranged blocks separated by a regular grid of 
fractures, wherein the transposed medium essentially 
exhibits a same fluid recovery during a process of 
capillary imbibition as the real medium. This being the 
case, the method consists in: 

- forming an image of the real medium in at least two 
dimensions in the form of an array of pixels, 

- determining for each pixel the minimum distance 
separating the pixel from the nearest fracture, 

- forming a distribution of numbers of pixels based on 
minimum distances to the fractures and determining the 
recovery function (R) of said set of blocks on the 
basis of this distribution and 

- determining dimensions (a, b) for equivalent regular 
blocks in the set on the basis of the recovery function 
<R) and the equivalent recovery function (Req) (using a 
procedure to identify said recovery functions, for 
example) . 
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With the method defined above based on a pixel 
representation of the medium, a number of different 
transfer fiinctions applied to Biny type of heterogeneous 
medium can be easily and rapidly confuted . 
5 The geometric method, for example, determines the 

dimensions of equivalent blocks, providing a very good 
match Of the imbibition behaviour of the real block or 
the distribution of real blocks, regardless of the 
shape (or shapes) of the blocks in question. Although 
10 simplified compared with the methods of the prior art, 
the oil recovery curve plotted for the equivalent 
section of the block has proved to be very close to 
that plotted for the real block. 

Other features and advantages of the method of the 
15 invention will become clearer from the following 
description of embodiments, given by way of example and 
not restrictive in any respect, and with reference to 
the appended drawings, of which; 

figure 1 illustrates a known procedure used to 
20 determine a medium having regular fractures equivalent 
to the real fractured medium, 

figure 2 illustrates the procedure of the 
invention which allows a medium having regular 
fractures equivalent to a real fractiured medium to be 
25 determined, 

figure 3 shows an example of adjacent pixels used 
to compute the value assigned to a pixel. 
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figure 4 is a histogram of a possible distribution 
of pixels in relation to distance from fractures, 

figure 5 is one possible variation in a normalised 
invaded area as a function of distance to fractures, 
5 figure 6 shows another exan^le of adjacent pixels 

in three different planes, Sk-i, sk and Sk-t-1 for a 3D 
computation of values attributed to a pixel, 

figure 7 shows how adjacent pixels can be extended 
in order to improve the computation of values 
10 attributed to a pixel and 

figure 8 shows, for the purpose of validation, how 
good a match is obtained for two oil recovery curves 
OR(t) determined using a real "comb- shaped" block on 
the one hand and an equivalent rectangular block on the 
15 other. 

A new simplified procedure which can be used to 
compute the dimensions of the block section equivalent 
to the "horizontal" section of a natural fractured 
medium is explained below. 

20 It should first of. all be pointed out that 

according to the "vertical" fracture hypothesis, i.e. 
perpendicular to the stratification planes, the matrix 
medium is continuous from one geological layer to 
another and the problem inherent in finding equivalent 

2 5 block dimensions becomes a two-dimensional one. 
Consequently, the problem faced in this instance is 
that of determining the equivalent square or 
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rectangular surface area of numerical matrix blocks for 
each layer or group of layers exhibiting similar 
fracture propeirties. 

Secondly, the eqpjivalence of a dual porosity model 
relative to a fractured model has to be established in 
terms of flow behaviour. The flows occurring within 
fractured oil reservoirs during field exploitation are 
essentially multi-phase and there are two main drive 
mechanisms as regards oil recovery from the matrix, 
capillary imbibition and gravity drainage. The effects 
of these two mechanisms are combined in the oil 
recovery methods which predominate as a strategy for 
developing many fractured reservoirs . Compositional 
mechanxBras such as diffusion represent another factor 
inherent in the methods used to recover gas. It is for 
this reason that the geometric method described below 
as a means of determining equivalent blocks is based on 
the concepts of multi-phase flow. 

An embodiment of the invention is described below 
with reference to figure 2, which consists in essence 
in matching the oil recovery function R(t) (of the 
natural fractured medium) obtained by the above- 
mentioned reference method with the known recovery 
function Req(t} of the tr€m.sposed medium in relation to 
a two-phase oil -water imbibition process (during a 
water-oil capillary imbibition drive process) , This 
matching is applied to each layer of the fractured 
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medium and then for sets of n layers. In this case, the 
recovery function R(t) obtained is the sum of the 
different functions Rn(t) of the n layers weighted by 
the corresponding thicknesses Hn. Since the fractures 
5 are vertical, only the hori2ontal dimensions of the 
equivalent block are determined. The task of matching 
the functions R(t) and Req(t) is therefore a two- 
dimensional one. 

1 Geometric fonmilation 

^0 Since the fractures are defined by the coordinates 

of their boundary points on a two-dimensional section 
XY of a layer, the imbibition process by which water is 
present in the fractures cind oil is present in the 
matrix blocks must be determined. It is assumed that 

15 the penetration of the matrix by the water is of the 
piston type. It is assumed that the function x-f(t) 
which links the advancing movement of the water front 
over time is the same for all the matrix blocks, 
irrespective of their shape, and for all the elementary 

20 blocks. Consequently, fitting the functions R(t) and 
Reqit) is equivalent to fitting the functions R (x) and 
Req{x). These functions physically define normalised 
zones penetrated by water as a function of the forward 
movement of the imbibition front through the fractured 

2 5 medium. 

In 2D, the analytical expression of Req<x) is as 
follows : 
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Rcg(x) = l 



-^(.-2.X.-2,) = 2{i.l>-^,^,e[o.^f.|); 




Re^(x) = l,x> 



where a and b are the dimensions of the equivalent 
5 rectangular or square block (a and b > 0) : 

The function R<x) does not have an analytical 
expression. It is computed by dividing up the section 



The section XY of the layer studied is considered 
as being an image in which each pixel represents a 
surface element. The pixels are regularly spaced apart 
at a pitch dx in the direction X and dy in the 

15 direction Y. The algorithm used is intended to 
determine for each pixel in this image the minimum 
distance which separates it from the nearest fracture. 

The image is translated by a table of two- 
dimensional real numbers: Pict [0 :nx+l . 0 ;ny+l] where nx 

2 0 and ny are the numbers of pixels in the image in the 
directions X and Y. In practice, the total number of 
pixels is in the order of one million, for example. The 
values of elements in the Pict table are the distances 
sought . 



XY of the layer studied in accordance with the 



algorithm defined below. 



10 



2 Algorithm used to connpute the function R (x) 
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Initialisation : All the pixels through which a 
fracture passes are at a distance zero from the closest 
fracture. The Pict table for these pixels is therefore 
initialised at value 0. This is done using a known 
algorithm (the Bresline algorithm, for example) to 
which are attributed the coordinates of pixels 
corresponding to the two ends of a fracture considered 
as a segment of a line and which initialises (at 0 in 
this specific case) the nearest pixels. The other 
elements of Pict are initialised at a value greater 
than the longest distance which exists between two 
pixels of the image. This value is nx.dx+ny.dy, for 
example . 

Coanputation : The required distance of a given 
pixel from the closest fracture is computed on the 
basis of the distance values already calculated for the 
adjacent pixels. A value is assigned to it which, if it 
proves to be lower than the value initially assigned to 
it, is the minimum of the values of adjacent pixels to 
which is added the distance of these pixels from the 
one being considered. 

This computation is performed in two successive 
phases. During the descending pass, the image is 
scanned line by line from top to bottom and from left 
to right (from Pict [1,1] to Pict [nx,ny)). The pixels 
taken into account therefore differ depending on 
whether a descending pass or an ascending pass is being 
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run. As can be seen from figure 3, the black and shaded 

pixels are those which will be taken into account 

respectively during the descending passes and ascending 

passes for the pixel Px. 

The oblique distance dxy being defined as : 
dxy = -Jdx^ ^dy^ , the algorithm is written 

where j=l to ny 
I where i-l to nx 

I I Pictti,j3= min Pict [i-1, j] +dx, : descending pass 

I I Pict [i-l J-i]+dxy, 

I t Pict Ei, j-i] +dy, 

I I Pict [i+l, j-l] +dxy, 

I I PictCi,j] 

I end of loop on i 

end of loop on j 

where j=ny to 1, 
i where i=lx to 1, 

I I Pict (i, j]»tninPict [i+1, j3+dx, : descending pass 

I I Pict[i+1, j+l]+dxy, 

I I PictCi, j+l)+dy, 

i I Pict[i-l, j+i]+dxy, 

I I Pict[i,jl 

I end of loop on i 

end of loop on j 
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Hlsbogram : From the table Pict calculated in this 
manner, a histogram can be plotted by classifying the 
non-zero values (those assigned to the pixels outside 
the fractures) by ascending order. 

The cumulated result of this histogram gives, for 
any distance delimiting two intervals of the histogram, 
the number of non-zero pixels whose value is less than 
this distance. When the method described is applied to 
a fractured porous medium where this distance 
corresponds to the advancing water front, the cumulated 
result of the histogram therefore indicates the area 
penetrated by water. The curve R{x) is obtained by 
dividing this cumulated result by the total number of 
non-zero pixels (to standardise it) . The number of 
intez-vals used on the abscissa of the histogram will 
correspond to the number of dividing points on the 
curve R(x) . It may be chosen as being equal to 500, for 
example . 

3 Finding the equivalent block dimensions 

At this stage, the function R(x) is known and it 
is now necessary to find the parameters [aj] 
(dimensions of the equivalent block which minimise the 
functional) : 



where N is the number of dividing points of R(x) and 
(Xi) are the abscissa of these dividing points. 
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Division along the ordinate of R(x) : 

In order to give the same weight to all the 
volumes of oil recovered during imbibition, the curve 
R(x) is re-divided at a constsmt pitch along the axis 
of the ordinates (figure 5) . The sequence (Xi) used by 
the functional is deduced from this division. 

Minimi cation of the ftmctional : 

Since a and b play symmetrical parts in the 
expression Req(a,b,x), the following functional is 
actually used: 

^(w,v) = ^(j?(x,)-Rc^(«,v,x,)) 
/-I 

J Re^(w,v,x) = uxjc + vxx^ 
where i , i , e 



Minimising this function amounts to finding the 
pair (w,v) for which 7'(w»v)=0. This is done by a Newton 
algorithm. 

Then, the pair (a,b) which is sought is derived 
from (m,v) . Three situations may arise: 

1) V > 0 meams that one of the values of the pair (a,b) 
is negative, which does not make sense physically 
speaking. v«0 is then incorporated in the expression of 
K^g{u,v,x). which means that the fractures are parallel. 
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The operation is repeated and the pair is computed 

as follows I 



2 

a = = 
u 

bin£ini 



-2 



2) The. case u +4v<0 does not make any physical sense 
either since it means that are not real. If we let 
u''+4v=o, this will impose a factor whereby the 
elementary block sought is square shaped <a«b) . After 
minimisation, the pair (a J) is computed as follows: 

- 4 
o = — 
w 

b = a 

3) For the other values of the pair (^,v), we have: 



1 



a = = 

V 

- -M V«* +4v 



4 Validation of the method of determining an 
equivalent block 

The geometric method outlined above for an 
imbibition transfer function was validated relative to 
a conventional and very costly reference method which 
requires a single-porosity multi -phase flow simulator 
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to divide the matrix blocks and fractures up in a 
manner allowing the recovery curves to be compared. 
Conventional two -phase flow simulations were run in 
order to validate the solutions provided by the 
5 geometric method described above. The validation may 
incorporate the following steps : 

a) computing the oil recovery function R(t) for the 
(real) geological section using the conventional method 
(reference solution) ; 

10 b) applying the geometric method to the real section, 

which gives a solution (a,b) ; 

c) using the conventional method again, computing the 
oil recovery function Req(t) on the equivalent block 
section having dimensions (a,b) determined beforehand 

15 and effecting a comparison with the oil recovery 
function of reference R(t), 

The geometric method gives equivalent block 
dimensions which ensdsle a very good match of the 
imbibition behaviour of the real block, regardless of 

20 the shape of the block being considered. The oil 
recovery curve computed on the equivalent block section 
is always very close to that calculated for the section 
of the real block, as illustrated in figure 8. 

The same validation was also successfully 

25 performed for a distribution of blocks exhibiting other 
different dimensions and shapes. 
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5 Otbex applications of tlie method 

Accxiracy of the coxbputed distances to fractures 

In the algorithm used to compute the distances of 
the pixels to the fractures on the basis of the 2D 
5 image, the accuracy of the computation can be improved 
by taking into accoxmt a larger number of pixels 
adjacent to the one being studied. Figure 7 shows the 
adjacent pixels taken into account when the influence 
of a pixel Px is extended to two lines and two columns 

10 before or after the pixel Px. Similarly, the black and 
shaded pixels are those which will be taken into 
account respectively during the descending passes and 
ascending passes, whilst those shown by a cross will be 
eliminated as they are redundant. 

15 To increase the accuracy of the computations still 

further, the rone of influence of the pixels could be 
extended still further <to 3 lines and 3 columns or 
more) . In practice, however, for the application 
presented above, such an extension does not impart any 

20 significant improvement to the final results. 

Extending the method to a three -dlmezisioiial object 
The algorithm givem above can be applied to a 
volume. In this case, each pixel represents an element 
of volume. The table Pict is replaced by a 3- 

2 5 dimensional table Pict3D [0 :nx+l . 0 :ny+l . 0 : nz+l] where 
nx, ny and nz are the numbers of pixels along X, Y and 
2. The adjacent pixels taken into account during the 
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descending and ascending passes in order to compute the 
horizontal plane number k at a pixel Px are illustrated 
in figure 6 . 



In the example that has been developed in relation 
to the study of two-phase imbibition (water-oil, for 
example) , the aim is to ^determine the sizes of the 
blocks in relation to the distance of points to the 

10 nearest fracture. The geometric method of the invention 
can also be used for other types of transfers between 
two contrasting media, such as transfers of heat 
between a well and a reservoir, for example. Above all, 
however, the "distance between pixels" function used in 

15 the algorithm above can be replaced by any function 
which links the points of an image. The value of this 
function between this pixel and its adjoining pixels 
taken into account for the purposes of the computation 
will have to be known for every pixel in the image. 

2 0 This function may express the transmissivity values 
between the grids of a reservoir, for example, the 
centres of which are the pixels of the image that will 
be used to calculate the drainage volume of a well 
bored through this reservoir. 

2 5 In one such situation, the two ascending and 

descending passes performed in the algorithm may not 
prove adequate for the purposes of finding a minimum 



5 



Extending to any function 
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value at any pixel of the image. The operation will 
therefore be repeated until the computed values no 
longer change . 

By us ing the notat ions presented above and 
5 assuming that the function F(i,j,k,l) will return the 
value of the function between the pixels (i,j) and 
(k,l),-the 2D algorithm becomes: 

change-true 
as long_as (change==true) 
10 change==f alse 

where j«=l to ny 
I where i«l to nx 
I I temp=Pict [i, j ] 

I I Pict [i, j] «min (Pict [i-1, j ] +F (i, j , i-1, j ) , :descending 
15 pass 

1 I Pict{i-1, j-l3 +F(i, j,i-l, j-1) 

I I Pictli, j-l]+F(i, j,i-l, j-1) 

I I PictCi+1, j-l3+F(i, j,i+l, j-1) 

I ( PictCiol 
2 0 II if ( Pi ct [i, j 3 otemp) changes true 

I end of loop on i 

end of loop on j 

where j =ny to 1,-1 
25 I where i=snx to 1,-1 
I I temp^Pict ti, j] 

I I Pict ti, j]=min{Pict [i+i, j3+F(i, j,i+i, j) , 
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: ascending pass 
Pict [i + 1, j+l]+P(io 
Pict ti, +Fti, j+1) 
Pict j+ll+F(i, j j+l> 
Pict ti, j] 
if (Pict [io , ] otemp change=true 
end of loop on i 
end of loop on j 
as long_as 
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CULXHS 

1. A method of simplifying the task of modelling 
an original/ heterogeneous, porous geological medium in 
the form of a transposed or equivalent medium such that 

5 the transposed medium is equivalent to the original 
medium relative to a given type of physical transfer 
function which is known in respect of the transposed 
Txiedium, the method consisting in: 

- forming an image of the geological medium in at least 
10 two dimensions (2D) in the form of an array of pixels 

and associating a specific initial value for said 
function to each pixel in the array, 

- determining step by step the value to be assigned to 
each pixel in said array for the physical transfer 

15 fimction, by reference to values of the function 
assigned to adjacent pixels in the image, and 

- determining a physical property of the transposed or 
equivalent medium by identifying the volume of the 
transfer function known for the transposed medium with 

20 the value of the transfer function determined step by 
step for the original medium. 

2. A method as claimed in claim 1, characterised 
in that the heterogeneous geological medium is crossed 
by an irregular network of fractures, all geometrically 

25 defined in blocks of irregular shapes and sizes. 

3. A method as claimed in one of claims 1 or 2, 
characterised in that said physical transfer function 
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represents a distance between different parts of the 
geological medium, 

4. A method as claimed in one of claims 1 or 2, 
characterised in that said physical transfer function 
represents transmissivities between different parts of 
the geological medium. 

S . A method as claimed in one of claims 1 or 2, 
characterised in that said physical transfer function 
represents thermal transfers between different parts of 
the geological medium such as transfers of heat between 
a reservoir and a well crossing through this reservoir. 

6. A method as claimed in one of claims 1 or 2, 
characterised in that said physical transfer function 
represents any transfer of mass or flow between 
different parts of the geological medium. 

7 . A method as claimed in claim 2 , used to 
determine, from an image of a real porous geological 
medium crossed by an irregular network of fractures, a 
transposed medium comprising a set of blocks arranged 
regularly and separated by a regular grid of fractures, 
said transposed medium essentially giving the same 
fluid recovery CReq) during a capillary imbibition 
process as the real medium, the method consisting in*. 

- forming an image of the real medium in at least two 
dimensions (2D) in the form of an array of pixels, 

- determining for each pixel the minimum distance 
separating the pixel from the closest fracture. 
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- forming a distribution of the number of pixels based 
on the minlTnum distance to the fractured medium and 
determining the recovery function (R) of said set of 
blocks on the basis of this distribution and 

- determining dimensions (a, b) for equivalent regular 
blocks in the trsuisposed medium on the basis of the 
recovery function (R) and the recovery function (Req> 
of the equivalent block. 

8 . A method as claimed in any preceding claim 
and substantially as hereinbefore described with 
reference to figures 2 to 8 of the accompanying 
drawings . 
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